Molecular dynamics simulation for heat transport in thin diamond nanowires 
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The phonon thermal conductivity in diamond nanowires (DNW) is studied by molecular dynamics 
simulation. It is found that the thermal conductivity in narrower DNW is lower and does not show 
obvious temperature dependence; a very small value (about 2.0 W/m/K) of thermal conductivity is 
observed in ultra-narrow DNW, which may be of potential applications in thermoelectric devices. 
These two phenomena are probably due to the dominant surface effect and phonon confinement 
effect in narrow DNW. Our simulation reveals a high anisotropy in the heat transport of DNW. 
Specifically, the thermal conductivity in DNW along [110] growth direction is about five times larger 
than that of [100] and [111] growth directions. The anisotropy is believed to root in the anisotropic 
group velocity for acoustic phonon modes in DNW along three different growth directions. 

PACS numbers: 65.80.-g, 62.23.Hj, 63.20.kg 
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I. INTRODUCTION 

Carbon-based materials have been attractive to both 
theoretical and experimental researchers for decades. 
From the bonding properties, these materials can be 
divided into two groups with different carbon-carbon 
bonds: sp2 bonding and sp3 bonding. The former in- 
cludes the three-dimensional graphite, zero-dimensional 
fullerenes such as Cgo, one-dimensional carbon nanotube, 
and two-dimensional graphene. These sp2 bonding mate- 
rials have been investigated one by one and have achieved 
some industrial applications. As well known, the sp3 
bonding diamond is the hardest material and of high 
thermal conductivity. Hsu et al. have successfully syn- 
thesized diamond nanowires (DNW)^, a new carbon- 
based material with sp3 bonding inner the system. DNW 
in the experiment are grown by a bottom-up technology, 
and the obtained samples are about 60-90 nm in diameter 
along the [110] growth direction. Babinec et. have fabri- 
cated 200-nm-diameter DNW by top-down method, and 
used the DNW samples as a single-photon source for the 
photonic and quantum information processing.^ After the 
gold rush for diamond, fuUerences, carbon nanotube, and 
graphene, the DNW is a competitive candidate of carbon- 
based materials for future novel applications. The DNW 
may inherit various advanced properties from diamond 
and gain new characters from its nanowire configuration 
of tunable surface to volume ratio (SVR). 

Prior to the experiment, in 2003, Barnard et al. have 
done a series of theoretical works to investigate the sta- 
bility and possible phase transition for the structure of 
DNWj^"— The ah intio calculations show that the struc- 
ture relaxation of DNW depends on both the surface 
morphology and the growth direction'^. The effect of 
the boron and nitrogen dopants on the structure stabil- 
ity of the DNW was also examined by the same group*. 
They also found that the energy band gap of DNW is re- 
duced considerably due to the effect of surface states^. 



In 2008, Tanskanen et.al carried out quantum chemi- 
cal calculations to study the the structural, electronic, 
and mechanical properties of DNW.— Since 2004, there 
is quite few theoretical works on the DNW, as it was 
a challenge to grow the DNW for experimentalists at 
that time. As now the DNW samples have been success- 
fully prepared, it is an urgent task to pursue theoretical 
studies on important properties of the DNW. The ther- 
mal transport is one of the very important phenomena 
in the quasi-one-dimensional nanowire structures, where 
the phonon spectrum modification and boundary scat- 



tering are important 
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Some possible practical appli- 



cations of the phonon confinement effect has also been 
foundi^ 

In this paper, we investigate the thermal conductiv- 
ity in the DNW by molecular dynamics simulation us- 
ing the Tersof&S and Brenner— empirical inter-atomic 
potentials. We show that the thermal conductivity of 
thicker DNW decreases with increasing temperature due 
to phonon phonon scattering (PPS), yet does not depend 
on temperature in narrower DNW where surface scat- 
tering is important. Our simulation discloses a highly 
anisotropic heat transport in DNW: thermal conductivity 
in DNW along [110] growth direction is about five times 
higher than that of [100] and [111] growth directions. By 
performing lattice dynamics analysis for DNW, we be- 
lieve that this anisotropy is originating from anisotropic 
group velocity of acoustic phonon modes in DNW along 
three different growth directions. 



II. STRUCTURE OF DNW AND SIMULATION 
DETAILS 

The DNW shown in Fig. [T] is denoted as 
DNW[110](60,10), indicating that DNW is growing 
along [110] direction with length L = 60 A and diameter 
Z? = 10 A. We adopt this notation in the following. The 
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FIG. 1: (Color online) Configuration of DNW[110](60,10). 
Two ends (blue online) are fixed. The high/low temperature 
controlled regions are on left and right. Other atoms (red 
online) are free during simulation. 



DNW is cut from the bulk diamond crystal by using a 
virtual cylinder with structure parameters (L, D)^, so 
that the cross section is always a circle and the diameter 
is well defined to be the diameter of the virtual cylinder. 
We applied free boundary condition for the surface of 
the nanowires, i.e atoms on the surface have free bonds. 
In thermal transport, these free bonds introduce some 
roughness to the surface of the nanowires and they 
can reduce thermal conductivity.— The macroscopic 
word 'roughness' is always borrowed to describe the 
microscopic structure of the nano-devices surface. The 
surface roughness can be introduced to manipulate the 
thermal conductivity through generating some vacancies 
on the surface which produces some free bonds on the 
surface. In this sense, the free boundary condition on 
the surface is actually equivalent to a rough surface. An 
experimental group has taken advantage of the rough 
surface to enhance the thermoelectric performance of 
silicon nanowires;^ Hence, the free boundary condition 
on the surface is of practical importance, so we prefer 
to use it in our simulation. In the molecular dynamics 
simulation, the free bonds of surface atoms are taken 
care naturally by the bond-order Tersoff potential. The 
local environment of surface atoms is different from inner 
atoms, and this difference is reflected by the bond-order 
parameters of the Tersoff potential. Both ends are fixed 
during MD simulation. The heat energy is pumped into 
the DNW through the temperature-controlled region 
on the left (black online) with thermal current Jl and 
will flow out of the DNW through the temperature- 
controlled region on the right (yellow online) with 
current J^. At steady state, energy conservation re- 
quires that Jl — — Jr. Using this relation, the thermal 
current flowing through the DNW can be obtained by 
J — (Jl — Jr)/2, with dJ — Jl + Jr. as an estimated 
of error. All other atoms (red online) vibrate freely. It 
is important to keep a free segment between fixed end 
and the temperature-controlled region to avoid energy 
accumulations and big boundary temperature jumpsJ^ 

The Nose-Hoove r^^i^^ thermostat is employed to 
maintain constant temperature in the left and right 
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FIG. 2: (Color online) Thermal quantities in 
DNW[110] (60,10) at 300K. (a). Temperature profile at 
three different simulation time. (b). Thermal current across 
the DNW. 



temperature-controlled regions. The force between car- 
bon atoms in the Newton equation is described by the 
commonly used Tersoff^'^ and Brenner^lA empirical inter- 
atomic potentials. The Newton equations are solved by 
velocity Verlet algorithm with an integration time step 
of 1.0 fs. A total simulation time is typically around 
2.0 ns and will be sufficiently extended to guarantee the 
achievement of steady state (if necessary) . Fig.[2]exhibits 
that the current difference dJ is almost zero at 1.0 ns and 
the thermal current has already reached a saturate value. 
There is no difference between temperature profiles at 2.0 
ns and 1.0 ns. These facts show that the simulation time 
is long enough. 



III. RESULTS AND DISCUSSION 

In the Nose-Hoover thermostat, the effective relax- 
ation time T determines the temperature fiuctuation for 
temperature-controlled regions and the strength of con- 
nection between heat bath and DNW. A small t yields 
small temperature fiuctuation due to strong interaction 
between heat bath and system, which facilitates energy 
injection into DNW. Thus the DNW can arrive at steady 
state faster and will save much simulation time. On 
the other side, large r corresponds to a weak heat bath 
with more fiuctuations as shown in Fig. EJb), resulting 
in longer simulation time. From Fig. [3l^a), all relaxation 
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FIG. 3: (Color online) Heat transport in DNW[110] (60,10) at 
300K with different relaxation time of heat bath. (a). Tem- 
perature profile, (b). Thermal current. 



time T=4.5, 1.4, and 0.4 ps accounts for almost the same 
temperature profile. We will use t = 1.4 ps in the fol- 
lowing, so that the DNW can reach steady state quickly 
and without introducing much thermal noise to system. 

Following the Fourier law, the thermal conductivity k 
is obtained as a ratio of current to temperature gradient: 
K = —J/{dT/dx). The temperature gradient is calcu- 
lated from linear fitting to the temperature profile in the 
free region between left and right heat baths as shown in 
Fig.dJ To drive the thermal current across the DNW, the 
left/right temperature-controlled regions in Fig.[T]are put 
in heat baths with /^j = (l±a)T. T is the average tem- 
perature of the DNW and a is the ratio of temperature 
difference. If a is large, the big temperature difference 
can drive a large thermal current and set up a big temper- 
ature gradient. As a result, the thermal conductivity will 
be almost unaffected provided that the Fourier law still 
valid. Fig. |4] confirm this result. Panels (a) and (b) show 
that the thermal current and the temperature gradient 
increase linearly with increasing a. The resulted thermal 
conductivity shown in panel (c) only display small vari- 
ation with a increasing from 0.05 to 0.9. We will keep 
a = 0.1 in following calculations. 

The thermal conductivity at room temperature in 
DNW[110] of different length and diameter are shown in 
Fig. m K of DNW[110] with length as 60 A and different 
diameters are displayed in panel (a). Only a small vari- 
ation in thermal conductivity is observed between thick 
DNW. The thermal conductivity decreases with decreas- 




FIG. 4: Heat transport in DNW[110] (60,10) at BOOK with 
different temperature ratio a. (a). The thermal current, (b). 
The temperature gradient, (c). Thermal conductivity. 



ing diameter, due to the enhancement of surface scatter- 
ing in narrower DNW. A very low thermal conductivity 
of about 2.0 W/m/K is obtained for the ultra-narrow 
DNW. The poor capability of heat transport may find 
some use in the thermoelectric applications where lower 
thermal conductivity is desired. However, we do not ob- 
serve the increase of thermal conductivity with decreas- 
ing diameter, which was discovered in silicon nanowires 
by Ponomareva et alr^. It should be note, however, from 
panel (a), the value of k almost does not change for the 
two DNW of largest diameters in the figure, but it is 
still difhcult to predict the value of thermal conductiv- 
ity for experimental DNW samples with about an order 
of magnitude larger in diameter. Panel (b) shows k in 
DNW[110] with same diameter 7 A. The thermal conduc- 
tivity increases with the increase of length below 40 A, 
However, no obvious change is found for k above 40 A, 
which indicates a diffusive thermal transport instead of 
ballistic transport. It is quite reasonable to observe dif- 
fusive thermal transport in DNW, since we have success- 
fully obtained the linear temperature profile from the di- 
rect MD simulation as shown in Fig. [2] (a) and Fig.jSja). 
As a result of surface effect, the value of thermal conduc- 
tivity in DNW is much smaller than that of the diamond 
(about 2000 W/m/K)^. Similar phenomenon has been 
found in silicon nanowiresiSi. 

Fig. ini shows the temperature dependence for thermal 
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FIG. 5: Thermal conductivity at 300 K in DNW[110] with 
different size (a). L = 60 A and (b). D = 7 A. 



conductivity in different DNW. For DNW[110](60,10), we 
have applied both Tersoff (circle, red online) and Bren- 
ner (square, red online) potentials to describe the in- 
teraction between carbon atoms. Both curves show a 
similar trend that k decreases with increasing tempera- 
ture in the whole temperature range. This result mani- 
fests the importance of the phonon-phonon scattering in 
DNW[110](60,10). The PPS happens inside the system 
and is stronger at higher temperature^^. The relaxation 
time of phonons corresponding to normal PPS at high 
temperature is inverse proportional to temperature T, 
leading to lower k at higher T. Compared with Ter- 
soff potential, the Brenner potential gives an obviously 
smaller value for n in whole temperature range. Both 
potentials are multi-body interaction and include non- 
linear interactions. They have been widely applied to 
describe the carbon-carbon interactions and the bond- 
ing process. It is not easy to compare the difference in 
the nonlinear interaction between these two potentials. 
However, we can extract some valuable information on 
the nonlinear properties of these two potentials from the 
thermal conductivity of these two potentials in Fig. [51 
With increasing T, two curves in the figure decline at 
almost the same speed, which implies that the nonlinear 
interactions in two potentials are more or less the same. 
Actually, the linear part of these two potentials has been 
compared, and the Brenner potential gives weaker inter- 
action than Tersoff potential for carbon systems. The 
frequency of optical phonon modes from Brenner poten- 
tial is overall smaller than that of the Tersoff potential. 



E 




[110](60,10) 
[1101(60,10) Brenner 
[110](60,7) 
[100](60,10) 
[111](60,10j^; 



900 



FIG. 6: (Color online) Thermal conductivity v.s temperature 
in different DNW with different types of interaction. 



This difference is responsible for the lower thermal con- 
ductivity from Brenner potential than Tersoff potential 
in Fig. [6] Experimental works are in need to judge which 
potential is more suitable to be applied in the study of 
thermal conductivity of DNW. We use Tersoff potential 
in all other calculations of this paper. 

It is quite interesting that the thermal conductiv- 
ity of DNW[110](60,7) in Fig. M does not depend on 
the temperature. This is different from the thicker 
DNW[110](60,10), where the thermal conductivity is de- 
termined by the PPS as discussed above. The thermal 
conductivity is not sensitive to temperature probably be- 
cause of the dominant boundary scattering and phonon 
confinement effects in thin DNW. In the thermal trans- 
port of nano-materials such as DNW, the surface effect 
is another important mechanism due to confinement ef- 
fect of the DNW. The PPS has dominant effect in bulk 
systems or nanowires with low SVR, while the surface ef- 
fect determines the thermal conductivity in ultra narrow 
DNW with high SVlU^i^. The surface scattering can be 
described in various forms. A more widely used formula 
for the relaxation rate of surface scattering can be found 
from Ziman's booki^ 
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where w^cr is the phonon velocity of mode (g, ct), and L is 
the sample size. Tbs is the corresponding life time and p 
is the spectacular parameter. The thermal conductivity 
can be obtained from the kinetic theory: 



(2) 



where V is the volumn and Cph(w) is the phonon heat 
capacity. The life time due to boundary scattering is tem- 
perature independent, so the temperature dependence for 
the thermal conductivity is the same as the heat capac- 
ity; i.e, K increases with increasing T. The tempera- 
ture dependence of k is determined by the competition 
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between the PPS and surface scattering. In thinner dia- 
mond nanowires, the surface scattering is also important; 
thus the interaction between these two mechanisms leads 
to much smooth temperature dependence for k and even 
results in the absence of \/T character. That is the ori- 
gin of temperature insensitive thermal conductivity in 
DNW[110](60,7). We have applied free boundary con- 
dition for the surface of the nanowire, so atoms on the 
surface have less than four neighboring carbon atoms. 
This particular surface configuration is taken care by the 
bond-order Tersoff potential which we have applied to de- 
scribe the interatomic interactions. During the molecular 
dynamics simulation, the Newton equations of all car- 
bon atoms are numerically solved, including the surface 
atoms. The movement of the surface atoms is different 
from the inner atoms, because of different local environ- 
ment in the Tersoff potential. Actually, there are some 
surface phonon modes on the surface, which contribute 
to the surface density of state as indicated by Fig. [71 As 
a result, the surface atoms introduce some roughness ef- 
fect to the thermal transport, which eventually leads to 
the surface scattering. The bond-order Tersoff potential 
includes the nonlinear interactions such as in the real 
space with u as the vibrational displacement. A Fourier 
transformation of this term leads to the three-phonon 
scattering process. So the PPS is also automatically in- 
cluded in the molecular dynamics simulation through the 
Tersoff potential. We can not distinguish surface scatter- 
ing and PPS in the molecular dynamics simulation in 
the real space, because these two scatterings are differ- 
entiated in the eigen mode space (the reciprocal space). 
The thermal conductivity depends strongly on the details 
of the surface configuration. From molecular dynamics 
simulation, different surface configuration will result in 
different local environment for surface atoms, which is 
reflected through the bond-order term of the Tersoff po- 
tential. This will directly affect the heat transport in 
the diamond nanowires. From kinetic theory, different 
surface configuration will result in different value of the 
spectacular parameter p in Eq. ([IJ . It also has straight- 
forward effect on the value of k. So there will be some 
difference if the surface boundary condition is not free as 
used in our work. 

Eq. ([T]) and ^ can give a more rigorous description for 
the surface scattering, yet there is an adjustable parame- 
ter p. This parameter need to be fitted to the experimen- 
tal resultsj^ As there is no experiment here, instead, we 
calculate the thermal conductivity k by the Debye model 
with formula: 

K = -^^CK,o-)%<T^ (3) 

^ a 

q,(j 

where Va is the volumn per atom. Cq^a and Vq,^ are 
the heat capacity and velocity of phonon mode {q,<j). 
q is the wave vector and a is the branch index for the 
phonon spectrum. In the quasi-one-dimcnsional diamond 
nanowire, the four acoustic phonon branches dominate 
the thermal conductivity, so only these four acoustic 
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FIG. 7: (Color online) Thermal conductivity v.s 1/diameter 
at 300 K in DNW[110] with length 60 A. The MD simulation 
results (black filled square) are compared with simple Debye 
thermal conductivity (blue line). 



phonon branches are considered. The velocities of these 
four acoustic branches are obtained through v = duj/dq 
from their spectrum. The mean free path I is taken 
to be the diameter of the nanowire. It corresponds to 
strong roughness of the surface, indicating the absorp- 
tion of phonon as soon as the phonon reaches the sur- 
face. This condition is reasonable for diamond nanowires 
of small diameter where strong boundary scattering has 
been found by the MD simulation. Fig. [7] shows the ther- 
mal conductivity from our MD simulation (black filled 
square) and this simple Debye model (blue line) . The di- 
amond nanowires studied here are DNW[110] with length 
60 A and different diameters. The four acoustic phonon 
velocities are 21.3, 10.8, 10.5, and 9.3 km/s. The Debye 
model is initially established in the bulk material; how- 
ever, it is quite interesting that this Debye thermal con- 
ductivity is close to the MD results for nanowires of small 
diameter, and k decreases with decreasing diameter, con- 
sidering its extraordinary simplicity. In the large diam- 
eter region, the Debye model overestimates the value of 
thermal conductivity, because the phonon-phonon scat- 
tering is completely ignored in this model. It should 
be note that the velocities for phonon modes near the 
boundary of Brillouin zone become smaller, which is an- 
other reason for the overestimation of k by Debye model. 

In Fig. [6l we compare the thermal conductivity of 
DNW in three different growth directions [100], [110], 
and [111]. These DNW are of the same size (60,10). A 
highly anisotropic thermal conductivity is observed: k in 
DNW[110] is about five times larger than the other two 
growth directions. To reveal the origin for this big differ- 
ence in thermal conductivity, we study phonon modes in 
DNW along three growth directions as shown in Fig. |S1 
In panels (a)-(c), the wave vector k is in the dimension 
of reciprocal lattice of the one-dimensional Brillouin zone 
space, where a = 2.51, 3.56, and 6.16 A is the lattice con- 
stant for DNW along [110], [100], and [111] growth direc- 
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FIG. 8: (Color online) Phonon modes in DNW with size 
(60,10) A in three growth directions, (a)-(c) are phonon 
dispersion for DNW in [110], [100], and [111] growth direc- 
tions, (d). The difference of phonon dispersion between three 
growth directions in low frequency region, (e). The DOS of 
phonon modes, (f). The DOS projected onto surface atoms. 



tions. We have shown all phonon modes in half of the full 
Brillouin zone for three DNW. The group velocity of the 
acoustic phonon modes can influence the thermal con- 
ductivity in three aspects. Firstly, a larger group veloc- 
ity directly leads to higher thermal conductivity, which 
is proportional to phonon group velocity and relaxation 
tim&i^. Secondly, larger group velocity will result in 
smaller possibility of Umklapp PPS (which can reduce 
thermal conductivity) for these acoustic phonon modes 
in low temperature region. This is because they have 
smaller wave vectors for the same frequency, thus the 
possibility of Umklapp PPS is reduced. This effect fur- 
ther enhances the thermal conductivity. The importance 
of Umklapp PPS on thermal conductivity was also found 
in graphene^i. Thirdly, with increasing wave vector fc, 
the dispersion curve for acoustic phonon mode of larger 
group velocity will increase and resonate with the disper- 
sion curve of optical phonon modes more quickly. The 
resonance between acoustic and optical phonon modes 
leads to strong PPS, thus reducing the thermal conduc- 
tivity in high temperature region. In one word, for sys- 
tems with acoustic phonon modes of larger group veloc- 
ity, the first two mechanisms will enhance the thermal 
conductivity, while the third mechanism reduce the ther- 
mal conductivity and leads to obvious temperature de- 
pendence of the thermal conductivity. The third mech- 
anism will compete with the first two mechanisms and 
it will be more important at high temperatures. For 



DNW, there are four acoustic phonon modes in the low 
frequency region: one longitudinal mode, two degener- 
ated transverse modes and a twisting mode. As the low- 
frequency phonon modes are important, we compare the 
four acoustic phonon dispersions in Fig. |6] (d). Obvi- 
ously, the group velocity of acoustic phonon modes are 
in the order of [110] > [100] > [111]. The largest group ve- 
locity of acoustic phonon modes in DNW[110] is due to 
its strongest Youngs moduli (22.8 GPa) that is about 
50% larger than the other two directions. The strong 
Youngs moduli may have relation to the redistribution 
of the elastic deformations during the relaxation of the 
surface atoms in the DNW, which happens in the acous- 
tically mismatched nanowires as discussed by Pokatilov 
et.al.^^-'^'^ As discussed above, the largest group veloc- 
ity in DNW[110] leads to highest thermal conductivity 
in DNW[110] and obvious temperature dependence for n 
shown in Fig. [51 Due to the same reason, the thermal con- 
ductivity is much smaller in DNW[100] and only depends 
on temperature above 500 K. The smallest group velocity 
in DNW[111] results in lowest thermal conductivity and 
temperature insensitive k in whole temperature range. 

Besides the low frequency acoustic phonon modes, we 
should examine the high frequency optical modes, which 
may play some role in the anisotropy of thermal con- 
ductivity even though these modes have limited contri- 
bution to K, due to their very small group velocity. We 
calculate the total density of state (DOS) for phonon 
modes in three growth directions as compared in panel 
(e). The major difference in DOS is found in high fre- 
quency region around 1600 cm~^. The DOS is obviously 
larger in [100] and [111] growth directions compared with 
[110] growth direction. At first glance, this should lead 
to higher thermal conductivity in these two directions. 
However, we further find that the surface modes are the 
major component of the high frequency phonon modes in 
DNW[100] and DNW[111]. This fact is disclosed by the 
surface DOS (SDOS), the projection of DOS onto surface 
atoms, as shown in panel (f). The SDOS is much larger in 
DNW[100] and DNW[111], resulting in very smaU contri- 
bution from these high frequency optical phonon modes 
to thermal conductivity in DNW[100] and DNW[111]. 



IV. CONCLUSIONS 

In conclusion, by carrying out molecular dynamics sim- 
ulation, we have investigated the thermal conductivity 
in the DNW with the widely used Tersoff and Bren- 
ner empirical inter-atomic potentials. We find that the 
thermal conductivity in thicker DNW is lower at higher 
temperatures, while k in narrower DNW is insensitive 
to temperature. The thermal conductivity shows de- 
creasing behavior with decreasing diameter of DNW and 
K « 2.0 W/m/K in ultra-narrow DNW which benefits its 
thermoelectric applications. We demonstrate that the 
thermal conductivity in DNWfllO] is about five times 
larger than that of DNW[100] and DNW[111] of the 
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same size. The obtained results are interpreted by PPS 
and the anisotropic group velocity of acoustic phonon 
modes. Diffusive transport is reached for DNW with 
lengths larger than 50 A and diameter of 7 A. 
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